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ABSTRACT 



We present optical I-band light curves of the gravitationally lensed double QSO 
B1600+434 from observations obtained at the Nordic Optical Telescope (NOT) between 
April 1998 and November 1999. The photometry has been performed by simultaneous 
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deconvolution of all the data frames, involving a numerical lens galaxy model. Four 
methods have been applied to determine the time delay between the two QSO com- 
ponents, giving a mean estimate of At = 51 + 4 days (95% confidence level). This is 
the fourth optical time delay ever measured. Adopting a = 0.3, A = Universe and 
using the mass model of Mailer et al. (2000), this time-delay estimate yields a Hubble 
parameter of Hq = 52^g 4 km s -1 Mpc -1 (95% confidence level) where the errors include 
time-delay as well as model uncertainties. There are time-dependent offsets between the 
two (appropriately shifted) light curves that indicate the presence of external variations 
due to microlensing. 

Subject headings: cosmology: observations — gravitational lensing: individual (B 1600+434) 
- distance scale — galaxies: spiral 



1. Introduction 

Intensive observational studies of gravitationally lensed QSOs have been conducted in the 
last few years with the aim of determining cosmological parameters, e.g., the Hubble constant, 
Hq, and to study the dark- matter distribution in lens galaxies. In particular, there has been a 
significant effort to measure time delays between lensed QSO components to derive Hq with the 
method described by Refsdal (1964). Following a more than decade- long monitoring of the double 
QSO 0957+561 by Schild and coworkers (Schild 1990; Schild & Thomson 1995), Schild & Thomson 
(1997) and Kundic et al. (1997) succeeded in pinning down its time delay in 1997. The same 
year saw an even more impressive accomplishment of Schechter et al. (1997) to determine multiple 
time delays in the 'triple QSO' PG1115+080 from a peak-to-peak variation of barely 0.15 mag. 
Following this demonstration of the feasibility of measuring time delays of multiply imaged QSOs 
an ongoing photometric monitoring was initiated at the Nordic Optical Telescope (NOT) in April 
1998. The program involves measuring time delays between lensed QSO components in as many 
systems as possible. Our main target during the first year of monitoring was the doubly imaged 
radio source QSO B1600+434 (Jackson et al. 1995) at redshift z = 1.59 (Fassnacht & Cohen 1998), 
which is gravitationally lensed by an edge-on late-type galaxy at z = 0.41 (Jaunsen & Hjorth 1997; 
Fassnacht k, Cohen 1998). The 1.39" angular separation between the two QSO images (labeled A 
and B in Fig. 1), together with the observed large flux variations in the system (Jaunsen & Hjorth 
1997), makes it well suited for time-delay measurements. Given the poor knowledge of the mass 
distribution in dark- matter halos of spiral galaxies, a time-delay measurement in B 1600+434 may 
not provide a firm determination of Hq. However, once Hq has been determined from other lensed 



1 Based on observations made with the Nordic Optical Telescope, operated on the island of La Palma jointly by 
Denmark, Finland, Iceland, Norway, and Sweden, in the Spanish Observatorio del Roque de los Muchachos of the 
Institute de Astrofisica de Canarias. 
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systems or using other methods, the time delay of B1600+434 can provide new constraints on the 
distribution of mass in the various components of spiral galaxies in general (Mailer et al. 2000). 

B1600+434 is very faint, with / ~ 22 for the faintest QSO component and / = 20.3 for the 
lensing galaxy. Furthermore, the B-component is substantially obscured by the lens (see Fig. 1) 
and photometry of the image is non-trivial. The data have therefore been analyzed with advanced 
deconvolution techniques both in order to model the light distribution of the lensing galaxy and to 
achieve accurate photometry of the QSO component blended with the lens nucleus. The time delay 
we present using the deconvolution technique is the fourth one measured at optical wavelengths, and 
it is the first of our monitoring program. A preliminary report was presented in Hjorth et al. (1999). 
The results presented here include more data points and supercede this report. Independently, a 
time delay has been measured at radio wavelengths with the VLA during the same observing season 
(Koopmans et al. 2000). 

2. Observations and Data Reduction 

Weekly observations of B1600+434 were carried out at the NOT from April 1998 to November 
1999, except for a short period (i.e., from December 1998 to February 1999) where the object 
was below the horizon at the NOT. Three different instruments were used: ALFOSC (Andaluci'a 
Faint Object Spectrograph), HiRAC (High Resolution Adaptive Camera) and the stand-by camera 
StanCam, equipped with detectors yielding pixel scales of / /188, / /107 and Of! 176 respectively. The 
/ band was chosen to minimize the effect of the extinction of the B component by the lensing 
galaxy. The total exposure times for each data point were adjusted according to the moon phase 
and typically ranged from 20 to 40 min, divided into three dithered exposures. The seeing usually 
varied from / /7 to l'/4, with 0?!9 being the most frequent value. A typical signal-to-noise ratio (S/N) 
of 100 was reached for the A component and 50 for the B component. Color terms were determined 
for the three detectors in order to match the light curves obtained with the three instruments. 

An automated pipeline employing routines in the IRAF/NOAO/CCDRED package has been 
developed in order to pre-process the CCD frames in an efficient and homogeneous way. Fringe- 
correction and cosmic-ray removal were performed on the individual frames before combination. 

3. Photometry 

The photometric data consist of one stacked frame per epoch. All light curves are calculated 
relative to 3 stars in the field that have calibrated / band magnitudes (Jaunsen h Hjorth 1997). 
Two of these stars, angularly close to the QSO, labeled SI and S2 by Jaunsen & Hjorth (1997), are 
used to construct the Point Spread Function (PSF). 

All the combined frames, one per epoch, are deconvolved with the MCS deconvolution algo- 
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rithm (Magain, Courbin & Sohy 1998). With this technique, many frames of a single object can be 
simultaneously deconvolved. This procedure combines the total S/N of all the frames obtained at 
different periods to determine the light distribution of the extended sources (galaxies) as well as the 
positions of the point sources (QSOs), since these parameters do not vary with time. The intensities 
of the point sources, however, are allowed to vary from one image to the other, hence producing 
the light curves. This technique is particularly well suited for the analysis of B1600+434 because 
of the light contamination of the faint B component by the spiral lensing galaxy. Simultaneous 
deconvolution of all the frames allows us to derive a high signal-to-noise numerical galaxy model 
(see Fig. 1). The derived galaxy is in good agreement with the images of B1600+434 obtained with 
the Hubble Space Telescope (HST) (Mailer et al. 2000). 

Maps of the residuals for each frame are used to check deconvolution results. These maps 
represent the x 2 fit m each pixel; i.e., the match between the deconvolved model image (reconvolved 
with the PSF) and the data. By inspecting these residuals we conclude that the photometry of 
both components was mainly limited by photon noise. 

In order to check for errors introduced by the PSF we deconvolve two fainter stars in the field 
that are not used in the construction of the PSF. Light curves for two of these stars, labeled S3 
and S4, are displayed in Fig. 2. Assuming that these stars are not variable, we use the standard 
deviation around the mean value as a measure of the photometric accuracy, giving la errors of 
0.026 mag and 0.039 mag for S3 and S4 respectively. For each point we subtract the photon noise 
in quadrature and attribute the residual error to the PSF. For the QSO images these estimated 
PSF errors are added in quadrature to the photon noise in order to model the total errors in the 
photometry. 

4. Light curves 

The /-band light curves (Fig. 2) contain 41 data points for each component. As predicted by 
lens theories, these light curves show that A is the leading component. A feature in A observed at 
JD 2451050 (September 1998) is repeated in B about 50 days later. In June 1999 (~ JD 2451300) 
there is an increase in the flux from both components. This increase is particularly strong in the 
B component and represents more than one magnitude in less than 100 days. After the intensity 
increase both A and B display several sharp short-term variations. Unfortunately our light curves 
are not sufficiently well-sampled during this period to recover the exact variations. This suggests 
that there are variations on time scales faster than the sampling frequency. 

The features in the light curves are sufficiently distinct that a rough eye-ball estimate directly 
yields an approximate time delay of about 50 days. Although 41 brightness measurements are a 
fairly small amount of data, application of standard techniques to measure time delays yield fairly 
robust results. We analyzed the light curves with the four methods described below. 
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5. Time delay measurements 

5.1. The SOLA method 

The method of subtractive, optimally-localized averages (SOLA) was originally developed for 
solving inverse problems in helioseismology (Pijpers & Thompson 1994). The basic idea of this 
method is to construct an optimal solution, taking into account measurement errors, of any linear 
inverse problem by taking linear combinations of the measured data. After having been applied 
successfully in helioseismology the SOLA method has found application in image reconstruction 
(Pijpers 1999), in the reverberation mapping of active galactic nuclei (Pijpers & Wanders 1994), 
and in the determination of time delays between lensed QSOs (Pijpers 1997). 

In the case of lensed quasar images, a transfer function which is a Dirac delta function is 
positioned at the time delay. With this method both the delay and the relative magnification of 
a pair of images are determined. Applied to B1600+434 a time delay of 55 days and a flux ratio 
of 0.79 is found. Also determined is a relative offset, for instance due to a contribution of the 
lensing galaxy or a foreground object which has not been subtracted, and a (linear) drift between 
the images which could occur due to microlensing events with long time scales. Problems can occur 
if there are higher order relative drifts such as short time-scale microlensing events which cannot 
easily be accounted for within this method which just uses determinations of low order moments 
of the transfer function. 

The errors are obtained by assuming that the photometric errors are uncorrelated and follow 
a Gaussian distribution. These propagated errors give a conservative error estimate of ±10 days 
on the time delay value, i.e., At = 55 ± 10 days. 

5.2. The minimum dispersion method 

For the combined data set generated from the A light curve and the shifted B light curve, 
we estimate the dispersion of the scatter around the unknown mean curve. The true time delay 
between the images should be manifested as a minimum in the dispersion spectrum (see Pelt et 
al. (1996) for definitions and notation). Both the simplest string length type statistic D\ and the 
smoothed dispersion spectrum D\ (smoothing parameter 5 = 10 days) give a clear global minimum 
at At pa 48 days. 

The precision of this time delay value is estimated using a bootstrap procedure. We first 
construct the combined light curve from A and the time-delay shifted B curve (At = 48). This 
yields a reference curve by median smoothing. From the reference curve we build bootstrap samples 
by adding randomized errors. For each sample (altogether 1000) we compute D| dispersion spectra 
and find corresponding dispersion minima. The scatter of the dispersion minima for the different 
bootstrap runs is significant (see Fig. 4), with a formal la error estimate of ±16 days. 
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5.3. Model fit method 

We model the data with an arbitrary continuous light curve with a fixed and constant sampling. 
For a given time delay value At, we use x 2 minimization to compare the model curve to the A 
curve and the appropriately time delay shifted and magnitude offset B curve. Additionally, a linear 
term in one of the curves can be included to model long term microlensing effects such as observed 
in QSO 0957+561 (Schild & Smith 1991). The minimization procedure is repeated for a range of 
time delay values (0 to 100) and the x 2 minimum is obtained with At = 49 days and Am = 0.66 
magnitudes (i.e., a flux ratio of 0.54). An additional small positive linear term for the A component 
further improves the x 2 value. 

A bootstrap method is used to estimate the errors. Two sets of 1000 curves are constructed 
with the same number of data points as in the data. One set is constructed with the same time 
sampling as that of our measured data, and the other with a random time sampling. Running the 
program on these simulated curves results in a la standard deviation of 2 days on the set of curves 
with the same sampling as our measured data, and of 7 days on the curves with random sampling. 
We could interpret this as an error of 2 days internal to the method since the sampling stayed the 
same for each simulated data set, and an error of 7 days on the measured time delay independently 
of the sampling of the data, but this needs to be investigated further. 

5.4. Iterative modeling 

Assuming that the additional time dependent offsets between the time-delay shifted curves are 
caused by microlensing, a fourth method based on iterative correction of the flux ratios is applied 
to the data. We split the light curves into several parts (bins) . In one case we chose three bins with 
lengths determined according to the apparent lengths of the offsets in a certain directions (positive 
or negative compared to the other curve). In three other cases we separate the curves into 2,3 and 
4 bins with an equal number of data points. For a range of time-delay values (40 - 65 days) we 
first fit a model (using the x 2 model fit method described above) independently to the separate 
bins. For each bin and time delay value, an additional constant or linear offset is determined to 
improve the fit. Finally, new time delay values are determined on these modified curves. The time 
delays determined in this way turn out to show very little sensitivity to the input values used, 
and to the different splitting of the curves. The results converge towards a mean value of 51 days 
with a standard deviation (la) of 2 days and magnitude offsets varying from 0.6 - 0.87 mag in the 
different bins. 
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5.5. Results 

The four time-delay estimates obtained from the different methods are consistent with one 
another (see Table 1). An average of the results gives a time delay estimate of 51 days and a flux 
ratio A/B = 1.50. The flux ratio at radio wavelengths has been measured to be 1.21 (Koopmans et 
al. 2000). If we assume that the radio flux is not affected by reddening, the B-component is reddened 
by 0.22 magnitudes, corresponding to a factor of 1.23. At a redshift of 0.41, this corresponds to 
Ay ~ 0.25 mag, essentially independent of the assumed reddening law. 

Except from the iterative modeling, the methods described above yield substantial errors in 
the time delay. However, the good agreement between the results in addition to the observed time 
dependent magnitude offsets in the time-delay shifted curves (Fig. 3) suggest that the large errors 
are mostly due to the presence of external variations, e.g., microlensing effects, rather than to the 
method used to determine the time delay. SOLA and the model fit method correct for slow, linear 
variations, but not higher order effects. Let us consider a case in which the brighter of the two 
lensed images undergoes many microlensing events during the time of monitoring. Such short term 
variations would be indistinguishable from measurement noise for these methods. They may fail to 
converge or produce a result that is biased at a level of the order of the errors. With the iterative 
method however, we also correct for higher order effects. The small errors from the iterative method 
compared to the other three methods indicate the presence of such higher order external variations 
in the light curves. 

If no assumptions are made on the short term offsets between the time delay shifted curves we 
must use the conservative error estimate of ±10 days (la) on the time delay value. However, if we 
assume that these external variations are higher order effects we can use the la error estimate of 
±2 days as found from the iterative method. 

6. Hq estimate 

Assuming that microlensing effects affect our I band light curves, a time delay of 51 ± 4 (95% 
confidence level) days has been estimated. Applying the galaxy models from Mailer et al. (2000), 
and assuming a 17 = 0.3 and A = Universe, our measured time delay is consistent with a Hubble 
parameter of Hq = 52^g 4 km s _1 Mpc -1 (95% confidence level) where both errors on the model 
and the time delay are included. 

In the model of B1600+434, Mailer et al. (2000) assume a constant M/L for the disk and 
the bulge, and a dark halo that has the same center and the same orientation as the disk. In 
addition they include the mass associated with the companion galaxy (see Fig 1) modeled as an 
isothermal sphere. The dark matter is modeled as a standard Pseudo Isothermal Elliptical Mass 
Distribution (Kassiola & Kovner 1993, PIEMD), and the exponential profiles of the disk and the 
bulge are modeled by a two-point PIEMD, also called the chameleon profile (Keeton & Kochanek 
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1998; Hjorth & Kneib 1999). Two types of solutions are found with these models, one with a large 
dark-matter halo core radius, and the other with a nearly singular dark halo. Applying only the 
solutions with a large dark halo core radius give H = 54^ km s' 1 Mpc- 1 . 

For comparison we also estimate Hq using a simple analytical generalized isothermal galaxy 
model as described by Witt et al. (2000). This method depends only on the time delay and the 
observed images positions. Using our measured time delay and the image positions from Mailer et 
al. (2000) measured on the HST image we obtain Hq = 45± 6 5 km s" 1 Mpc" 1 (Q = 0.3 and A = 0) 
where errors come only from the time delay and the image positions. This is only an indicative 
estimate of Hq using a very simplified model. Possible shear from ellipticity in the lens galaxy or 
mass contribution from the companion galaxy is not taken into account. 

Finally, as pointed out by Romanowsky & Kochanek (1999) and Witt et al. (2000), the time 
delay depends on the density profile of the lens galaxy. Hence if the real density profile of the spiral 
galaxy is different from our models the errors in Hq will increase. 

7. Discussion 

The main goal of our optical monitoring program at the NOT is to measure time delays for as 
many lensed QSOs as possible in order to better constrain the mass distribution in lens galaxies. A 
robust time delay value of At = 51 ± 10 (1 a) days and a flux ratio of 0.69 have been measured from 
/ band light curves of our first target B 1600+434. Assuming that our lightcurves are affected by 
external variations the time delay value is constrained to 51 ± 4 days (95% confidence level). This 
value is consistent with a Hubble parameter Hq = 52^g 4 km s _1 Mpc -1 (f2 = 0.3, A = 0) using the 
models of Mailer et al. (2000). We recall that possible systematic errors due to uncertainties in the 
density profiles and the degeneracy between the relative contribution of disk, bulge and halo to the 
mass in the spiral galaxy may increase the uncertainties in the estimated Hq. 

As is evident from Fig. 3, and as measured by the various methods to determine the time 
delay, there are slow additional time dependent magnitude offsets between the curves. These 
variations are ~0.2 mag on time scales of a few months. Since the QSO images of B1600+434 
pass through the lens galaxy with high optical depths for microlensing, the offsets between the two 
time-delay shifted curves (Fig. 3) are likely to be due to microlensing of one of the components. 
This interpretation is supported by the detection of microlensing in light curves obtained at radio 
wavelengths (Koopmans et al. 2000), and implies that a significant fraction of the mass in the 
dark halo consists of massive compact objects. Simulations of microlensing light curves must 
be carried out in order to determine the lens masses and source sizes that could reproduce the 
observed variations in our optical light curves (see Wambsganss & Paczyhski (1991) and Schmidt 
& Wambsganss (1998)). Such a microlensing analysis will be published in a separate paper. 

Microlensing could also be partly responsible for the sharp event at ~JD2451350 in the light 
curve of the B component. Such high magnification and short duration events may occur in cases 
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of microlensing by stars in random motion whereas the slow variations are typical for microlensing 
events by the stars with velocities following the bulk motion of the galaxy (e.g., Wambsganss & 
Kundic 1995). We note however that much of this event must be due to intrinsic variations in the 
QSO since a significant peak is detected for both components. 

Although the temporal sampling of our curves does not allow to fully disentangle high frequency 
microlensing and intrinsic variations, it is sufficient to follow lower order variations, yielding robust 
time delay estimates fairly independently of the statistical method used. The time delay estimates 
measured with the four different methods are in agreement with each other. Furthermore, the time 
delay value of 47tg 2 days recently estimated from radio measurements Koopmans et al. (2000) 
is consistent with our optical measurement. There are thus two independent data sets, one at 
optical wavelengths (J band) and one at radio wavelengths, hence affected in very different ways 
by microlensing, yielding the same time delay value within the measurement errors. 

The remarkably strong variations in this system makes it interesting to study in more detail. 
The observed microlensing effects may provide important constraints on the MACHO masses in 
the lens galaxy. Moreover, once Hq is known independently from other lenses, or other methods, an 
accurate time delay can be used to constrain the mass distribution between the halo and the bulge 
of the spiral lens galaxy. In particular, a well determined time delay and flux ratio can contribute 
to determine the maximum allowed disk mass in the lens galaxy. 

We thank the NOT Director Vilppu Piirola for granting us observing time for this project on 
a flexible basis. We are especially grateful to the dozens of visiting observers at NOT who have 
contributed to this project by performing the scheduled observations. This project was conceived 
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We thank Rolf Stabell and Sjur Refsdal for inviting us there and for their kind hospitality. We 
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the referee. The project was supported in part by the Danish Natural Science Research Council 
(SNF). IB was supported in part by contract ARC94/99-178 "Action de Recherche Concertee 
de la Communaute Frangaise (Belgium)" and Pole d'Attraction Interuniversitaire, P4/05 (SSTC, 
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Fig. 1. — Top: Stacked /-band images of 12 x 12 arcsec around B1600+434 with a total of ~ 3.5 
hours of exposure and a seeing FWHM = l'/13. Middle: The image (FWHM = 0'.'38) as obtained 
from the simultaneous deconvolution of 33 frames. The main lensing spiral galaxy can be seen 
between the two QSO components A and B. A neighbor galaxy is seen to the south east of the 
system. North is up and East is to the left. Bottom: The HST H band image. 

Fig. 2. — /-band light curves for B1600+434 (A and B component) and two reference stars in 
the field. The plotted value for B is /(mag)+1.5, and for S4 /(mag)+0.4. The magnitudes are 
calculated relative to calibrated stars in the field. The error bars represent photon noise and PSF 
errors measured as described in the text (§3). 

Fig. 3. — Combined lightcurve from both components of B 1600+434. The curve from the B 
component is shifted forward in time by 51 days and scaled with —0.69 mag. The magnitudes are 
calculated relative to calibrated stars in the field. The error bars represent photon noise and PSF 
errors measured as described in the text (§ 3). 

Fig. 4. — Distribution of the time delays for 1000 bootstrap runs with the minimum dispersion 
method (see § 5.2). 



Table 1. Estimated time delays and /-band magnitude differences for B1600+434 calculated 
with the four methods described in the text (§ 4). 



SOLA 
Minimum dispersion 

x 2 fit 

Iterative fit 



At (days) Am (mag) 



55±10 0.72±0.007 

48±16 

49±7 0.66±0.01 

51±4 0.6-0.87 



- 12 - 



B 1600+434 (NOT I) 

1 arc 


A 

sec 


Deconvolved 

• 

1 arcsec 


HST(H) 
■ 

m 

W 

1 arc 


« 

sec 



-13- 



20 



21 



cd 

£ 



23 



-1 




1 i 

1 A 


i i 




I 




i 




1 


- A A $ 




A A A A^ A A A A 
A ^ 

A 


AA 


A 


A 


A^A^A 


A AA. 


AA 
AA 


S4 . 


- 




• 


•• 


• 


• 


v— 

m 


• 




S3 _ 

- 


- 










03 








- 


















— 




CD 




03 




in 


03 m 
CD 




— 


03 














ffi 




A _ 

- 


- 03 ™03 
03 

- 


m 


03 

03 

03,- 


03 
03 










0# 
03 






% 








• 

* 

■ 




















i 


i 




B _ 














i 


i i 














i 














1 

i 














lit 


i 


1 4 


i 








«I 






i 




i i 




1 


i 


i 







200 400 600 

JD - 2450900 



-14- 



• A 

20 B-0.69 

shifted 51 days 



20.5 



21 - - 



21.5 



22 







200 



JD 



400 

2450900 



600 



-15- 




-40 -30 -20 -10 



10 20 30 40 50 60 70 80 90 
Delay (days) 



